#!/usr/bin/env Rscript

rm(list=ls())
library(FactorHet)
library(dplyr)
library(glue)


repdata <- readRDS('packaged_data.RDS')

fac<-c("Ed", "Gender", "Country", "Reason", 
       "Job", "Exp", "Plans", "Trips", "Lang")

hhy_fmla <- as.formula(glue::glue('Chosen_Immigrant ~ {paste(fac, collapse = " + ")}'))
int_hhy <- update.formula(hhy_fmla, '. ~ . + (Job + Ed + Reason) * Country +Job*Ed')

out <- list()
for (K in 1:4){
  
  set.seed(78712)

  #Run model
  est <- FactorHet_mbo(formula = int_hhy, K = K,
    design = repdata, 
    moderator =  ~ scale_hisp_prej_flip + party_ID + 
      ppEducat + census_div + ppEthm,
    group = ~ CaseID, task = ~ contest_no, choice_order= ~ choice_id)
  out[[K]] <- est
}

saveRDS(c(out, list(info = sessionInfo())), 'out_BIC.RDS')


out <- list()
for (K in 1:4){

  set.seed(78712)

  #Run model
  est <- FactorHet_mbo(formula = int_hhy, K = K,
    design = repdata, 
    moderator =  ~ scale_hisp_prej_flip + party_ID + 
      ppEducat + census_div + ppEthm,
    group = ~ CaseID, task = ~ contest_no, choice_order= ~ choice_id, 
    mbo_control = FactorHet_mbo_control(criterion = 'AIC')
  )
  out[[K]] <- est
}

saveRDS(c(out, list(info = sessionInfo())), 'out_AIC.RDS')
